
[FImgname,Pname] = uigetfile('*.txt','Choose data file');

s=strcat(Pname,FImgname);
f1=fopen (s,'r');



m=12;                                % number of curves

n=1216;                              % number of time frames in a file

data=zeros(m,n);                    % data (second column)
timer=zeros(n,1);                   % timer (first coloumn)


% reading data from file

plchar=zeros(1,m);
caarea=zeros(1,m);
mearea=zeros(1,m);

for j=1:n
    
    timer (j)=fscanf (f1,'%g',1);
    
   if j==1
             
    for i=1:m
        
        caarea(1,i)=fscanf (f1,'%g',1);
        cadat (i,j)=fscanf (f1,'%g',1);
        mearea(1,i)=fscanf (f1,'%g',1);
        mepdat (i,j)=fscanf (f1,'%g',1);
         plchar(1,i)=fscanf (f1,'%g',1);
    end;

   else
      for i=1:m
        caarea(1,i)=fscanf (f1,'%g',1);
        cadat (i,j)=fscanf (f1,'%g',1);
        mearea(1,i)=fscanf (f1,'%g',1);
        mepdat (i,j)=fscanf (f1,'%g',1);
             
        end; 

    end;
end;




for ncur=1:m                                         % !!!! number of platelet to analyze 



me=mepdat(ncur,:);

f=figure;
figure(f);
a=axes;
axes(a);
plot (me);
set (a,'nextplot','add');


prompt = 'How many release acts do you see here? ';
str = input(prompt,'s');

numreleases=str2num(str);
real_secretion_frame=0;


if (numreleases>0)

'Choose 4 points for the first release start and end '

[x,y]=ginput(4);

tstart=timer(round(x(2)))
tend=timer(round(x(3)));

mebefore=0;
j=0;
for i=round(x(1)):round(x(2))
    mebefore=mebefore+me(1,i);
    j=j+1;
end;
mebefore=mebefore/j;

x1=round(x(1)):round(x(2));
y1=x1*0;
y1=y1+mebefore;

h=plot (x1,y1);
set (h,'color',[1 0 0]);
set (h,'linewidth',2);



meafter=0;
j=0;
for i=round(x(3)):round(x(4))
    meafter=meafter+me(1,i);
    j=j+1;
end;
meafter=meafter/j;

x1=round(x(3)):round(x(4));
y1=x1*0;
y1=y1+meafter;

h=plot (x1,y1);
set (h,'color',[1 0 0]);
set (h,'linewidth',2);


real_secretion_frame=round(x(2));          % !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!  number of frame when release detected
real_secretion_endframe=round(x(3));

end;


scur='mepacrine';
scur=strcat(scur,'_',num2str(ncur),'.tif');
saveas (f,scur);





ca=cadat(ncur,:);

% zz=[0:max(ca)/100:max(ca)];
% 
% yy=histc(ca,zz);
% 
% f=figure;
% figure(f);
% 
% plot (zz,yy);
% 
% chk=0;
% 
% jsp=0;
% for i=2:length(zz)-1
%     
%     if ((yy(1,i)<yy(1,i+1)) && (yy(1,i)<yy(1,i-1)) && (chk==0) )
% 
%             jsp=i;
%             chk=1;
%     end;
% end;
% 
% 'min in hist'
% zz(1,jsp)
% 
% ca2=ca;
% for i=1:length(ca)
%     if ca(1,i)>zz(1,jsp)
%         ca2(1,i)=0;
%     end;
% end;
% 
% f=figure;
% 
% figure (f);
% 
% plot (ca2);
%     
% 
% 



% 
% 
% 
% 
% 
% 
checker=0;

while (checker==0)
    
f=figure;
figure(f);
a=axes;
axes(a);
h=plot (ca);
set (a,'ylim',[0 max(ca)] );
set (a,'nextplot','add');

prompt = 'How many baseline dots would you like to enter? ';
str = input(prompt,'s');

nn=str2num(str);

'Choose the points to build the baseline '

[x,y]=ginput(nn);

x;
y;


c1=polyfit (x,y,1);

naklon=c1(1,1);

y0=c1(1,2);

% !!!!!!!!!!!!!            %  linear  slope  

baseline=zeros(1,n);
for i=1:n
    baseline (1,i)=y0+(naklon*i);
end;

h2=plot (baseline);
set (h2,'color',[1 0 0]);


scur='calcium_baseline';
scur=strcat(scur,'_',num2str(ncur),'.tif');
saveas (f,scur);



c3=polyfit (x,y,3);

x3=1:n;

y3=(c3(1,1)*x3.^3) + (c3(1,2)*x3.^2) + (c3(1,3)*x3) + c3(1,4);

h3=plot (x3,y3);

set (h3,'color',[0 1 0]);





prompt = 'Good? y/n [y]: ';
str = input(prompt,'s');
if isempty(str)
    str = 'Y';
end

if ((str=='Y') || (str=='y'))
    
    checker=1;
else
    checker=0;
    
end;

end
    


dif=ca-baseline;



checker=0;

while (checker==0)

f=figure;
figure(f);

a=axes;
axes(a);
h=plot (dif);

set (a,'ylim',[0 max(dif) + max(dif)/10] );
set (a,'nextplot','add');


'Select critical Ca level to define the Calcium peaks'

[x,y]=ginput(1);


level1=y                                        %!!!!!!!!!!!             only values above this level are considered during peak search



levelmas=ones(1,n)*level1;

plot (levelmas);

peaks=zeros(1,n);

% simple way to search for a peak (local maximum condition for 5 peaks )

for i=3:n-2
    if  ((dif(1,i)>level1) && (dif(1,i-1)<dif(1,i)) && (dif(1,i+1)<dif(1,i)) && (dif(1,i+2)<dif(1,i)) && (dif(1,i-2)<dif(1,i)))
        
        peaks(1,i)=dif(1,i);
    end;
end;

masi=1:n;

scatter (masi,peaks);     % to plot peaks onto primary data 


prompt = 'Good? y/n [y]: ';
str = input(prompt,'s');
if isempty(str)
    str = 'Y';
end

if ((str=='Y') || (str=='y'))
    
    checker=1;
else
    checker=0;
    
end;



end;






sumca=zeros(1,n);         % calcium integral 


for i=2:n
    
    
    sumca(1,i)=sumca(1,i-1)+dif(1,i);
    
end;
    







% calculation of calcium integration profiles (for various sizes of integrating
% window 


okno=3;

sumca3=zeros(1,n);

for i=okno:n
    
    for j=i-okno+1:i
    sumca3(1,i)=sumca3(1,i)+dif(1,j);
    end;
    
end;

okno=5;

sumca5=zeros(1,n);
for i=okno:n
    
    for j=i-okno+1:i
    sumca5(1,i)=sumca5(1,i)+dif(1,j);
    end;
    
end;

okno=10;
sumca10=zeros(1,n);

for i=okno:n
    
    for j=i-okno+1:i
    sumca10(1,i)=sumca10(1,i)+dif(1,j);
    end;  
end;


okno=20;
sumca20=zeros(1,n);
for i=okno:n
    
    for j=i-okno+1:i
    sumca20(1,i)=sumca20(1,i)+dif(1,j);
    end;
    
end;


% f=figure;
% figure (f);
% 
% plot (sumca);


scur='spikes_data_new';
scur=strcat(scur,'_',num2str(ncur),'.txt');
f2=fopen(scur,'w');        % writing data to the file : number of peak, time of peak, amplitude of peak, current value of calcium integral

k=0;

% writing level of baseline 
    fprintf (f2,'%g    ',level1);
    % writing the absolute value of intensity which is further considered
    % as zero (initial value of ca intensity of resting platelet )
    
    fprintf (f2,'%g     ',y0);

    
         % writing  number of granule release acts 
          fprintf (f2,'%g    ',numreleases);
          
           
            % writing  platelet type (spreaded or non-spreaded)
          fprintf (f2,'%g    ',plchar(1,ncur));
          
           % writing  platelet area (mepacrine) )
          fprintf (f2,'%g    ',mearea(1,ncur));
          
            % writing  platelet area (calcium)
          fprintf (f2,'%g    \n',caarea(1,ncur));
          



for i=1:n
    if (peaks(1,i)>0)
        k=k+1;
        fprintf (f2,'%i %g %g %g \n',k,timer(i),peaks(1,i),sumca(1,i) );
    end;
end;

maspeak=zeros(1,k);
maspeaktime=zeros(1,k);
j=0;
for i=1:n
     if (peaks(1,i)>0)
         j=j+1;
        maspeak(1,j)=peaks(1,i);
        maspeaktime(1,j)=timer(i);
     end;
end;

maspeak;
avpeak=sum(maspeak)/k;
avpeak;


fclose (f2);




% writing integrated calcium profiles for windows of 3, 5, 10 and 20 points





sintegrat=strcat ('integrated_calcium_profiles_full_',num2str(ncur),'.txt');
f2=fopen (sintegrat,'w');

for i=1:n
   
        fprintf (f2,'%g %g %g %g %g \n',dif(1,i),sumca3(1,i),sumca5(1,i),sumca10(1,i), sumca20(1,i) );
   
end;

fclose (f2);


actframe=0;
for i=1:n
    if (timer(i)==0)
        actframe=i;
    end;
end;


%%%%        writing data for most recent peak (before granule release) to the file 

if (real_secretion_frame>0)


z=real_secretion_frame;
check=0;
jpeak=0;

for j=z:-1:1
    if ((peaks(1,j)>0)  && (check==0))
        check=1;
        jpeak=j;
    end;
end;

j=0;

for i=1:jpeak
    if (peaks(1,i)>0)
        j=j+1;
    end;
end;

peaknumber=j;

j=0;
sumpeaks=0;
for i=1:jpeak
    if ((peaks(1,i)>0) && (timer(i)>0))
        j=j+1;
        sumpeaks=sumpeaks+peaks(1,i);
    end;
end;

peaknumber_after_act=j;
avpeakaft=0;
    if (j>0)
        avpeakaft=sumpeaks/peaknumber_after_act;
    end;

    
j=0;
sumpeaks=0;
for i=1:jpeak
    if ((peaks(1,i)>0) && (timer(i)<0))
        j=j+1;
        sumpeaks=sumpeaks+peaks(1,i);
    end;
end;

avpeakbef=0;
if (j>0)
    avpeakbef=sumpeaks/j;
end;

peaks_before=j;

time_after_act=timer(real_secretion_frame);
time_full=timer(real_secretion_frame)-timer(1);

freq_before=peaks_before/abs(timer(1));
freq_after=peaknumber_after_act/time_after_act;




smain=strcat(Pname,FImgname);


smain2='a';
for i=1:length(smain)-4
    smain2(i)=smain(i);
end;

smain3=strcat(smain2,'_',num2str(ncur),'.txt');

f2=fopen (smain3,'w');


if (jpeak==0)
    
    fprintf (f2,'no spike before release detected! ' );
else
    
  
    % writing level of baseline 
    fprintf (f2,'%g    ',level1);
    
    % writing the absolute value of intensity which is further considered
    % as zero (initial value of ca intensity of resting platelet )
    
    fprintf (f2,'%g     ',y0);
 
    % writing  average level of Ca spikes
      fprintf (f2,'%g    ',avpeak);
    
    % writing  time of first release start
      fprintf (f2,'%g    ',tstart);
      
      
     % writing  duration of first release
      fprintf (f2,'%g    ',tend-tstart);
      
      % writing  amplitude of first release
      fprintf (f2,'%g    ',mebefore-meafter);
      
       % writing  number of the recent peak (from begining)
      fprintf (f2,'%g    ',peaknumber);
      
       % writing  number of the recent peak (from activation)
      fprintf (f2,'%g    ',peaknumber_after_act);
      
        % writing  cal frequency before activation (per minute)
      fprintf (f2,'%g    ',freq_before*60);
      
        % writing  cal frequency after activation (per minute)
      fprintf (f2,'%g    ',freq_after*60);
      
          % writing  av cal peak value before activation 
      fprintf (f2,'%g    ',avpeakbef);
      
       % writing  av cal peak value after activation (before release)
      fprintf (f2,'%g    ',avpeakaft);
      
        % writing  cal integral before release (full)
          fprintf (f2,'%g    ',sumca(1,real_secretion_frame));
          
          % writing  cal integral before release (after activation)
          fprintf (f2,'%g    ',sumca(1,real_secretion_frame)-sumca(1,actframe));
          
          
          % writing  cal integral at release ( averaged by 3 frames before release)
          fprintf (f2,'%g    ',sumca3 (1,real_secretion_frame) );
          
          
          % writing  cal integral at release ( averaged by 5 frames before release)
          fprintf (f2,'%g    ',sumca5 (1,real_secretion_frame) );
          
          
          
          % writing  cal integral at release ( averaged by 10 frames before release)
          fprintf (f2,'%g    ',sumca10 (1,real_secretion_frame) );
          
          
          % writing  cal integral at release ( averaged by 20 frames before release)
          fprintf (f2,'%g    ',sumca20 (1,real_secretion_frame) );
          
          
              % writing  amplitude of ca peak preceding release
          fprintf (f2,'%g    ',peaks(1,jpeak));
          
                % writing  time of ca peak preceding release
          fprintf (f2,'%g    ',timer(jpeak));
          
        % writing  delay time between first release and preceding ca peak
          fprintf (f2,'%g    ',timer(real_secretion_frame)-timer(jpeak));
          
           % writing  number of granule release acts 
          fprintf (f2,'%g    ',numreleases);
          
           
            % writing  platelet type (spreaded or non-spreaded)
          fprintf (f2,'%g    ',plchar(1,ncur));
          
            % writing  platelet area (mepacrine) )
          fprintf (f2,'%g    ',mearea(1,ncur));
          
            % writing  platelet area (calcium)
          fprintf (f2,'%g    ',caarea(1,ncur));
          
          
          
          
      
     
     fclose (f2);
end;


xpeak=jpeak;
ypeak=peaks(1,jpeak);

h3=scatter(xpeak,ypeak);
set (h3,'MarkerFaceColor','b');

y2=[0 1.5*ypeak];
x2=[real_secretion_frame real_secretion_frame];

h4=plot (x2,y2);

set (h4,'color',[1 0 0]);
set (h4,'linewidth',2);


scur='calcium_and_release';
scur=strcat(scur,'_',num2str(ncur),'.tif');
saveas (f,scur);



else 
    
smain=strcat(Pname,FImgname);

smain2='a';
for i=1:length(smain)-4
    smain2(i)=smain(i);
end;

smain3=strcat(smain2,'_',num2str(ncur),'.txt');

f2=fopen (smain3,'w');



    % writing level of baseline 
    fprintf (f2,'%g    ',level1);
    
    % writing the absolute value of intensity which is further considered
    % as zero (initial value of ca intensity of resting platelet )
    
    fprintf (f2,'%g     ',y0);
 
    % writing  average level of Ca spikes
      fprintf (f2,'%g    ',avpeak);
    
         % writing  number of granule release acts 
          fprintf (f2,'%g    ',numreleases);
          
           
            % writing  platelet type (spreaded or non-spreaded)
          fprintf (f2,'%g    ',plchar(1,ncur));
          
           % writing  platelet area (mepacrine) )
          fprintf (f2,'%g    ',mearea(1,ncur));
          
            % writing  platelet area (calcium)
          fprintf (f2,'%g    ',caarea(1,ncur));
          
          
     
     fclose (f2);
     
     end;


if (numreleases>0)
    
    
    
%%%%%%          GENERATING RANDOM SECRETION TIME and writing the data for
%%%%%%          most recent peak (for stasticial analysis)






tries=100;



scur='spike_data_random_secretion';
scur=strcat(scur,'_',num2str(ncur),'.txt');

f2=fopen(scur,'w');


for i=1:tries
    
    
    
cnt=0;
zcur=n;

while cnt<numreleases
    
z=round((rand(1)*(n-1)))+1;
if (timer (z)>0)
    cnt=cnt+1;
    if (z<zcur)
        zcur=z;
    end;
end;
end;

z=zcur;

check=0;
jpeak=0;



for j=z:-1:1
    if ((peaks(1,j)>0)  && (check==0))
        check=1;
        jpeak=j;
    end;
end;

jpeak





if (jpeak==0)
    
    fprintf (f2,' no spike before release detected! \n' );
else
    
    %   writing pseudo-realease time, recent peak time, delay, peak value, current value (in the moment of release) of full calcium integral, summ of previous 3 ca values, 5 values, 10 values and 20 values
    
   fprintf (f2,' %g %g %g %g %g %g %g %g %g \n', timer(z),timer (jpeak), timer(z)-timer(jpeak), peaks(1,jpeak), sumca(1,z), sumca3(1,z),sumca5(1,z),sumca10(1,z), sumca20(1,z) )
    
    
     
end;
    
    
end;

fclose (f2);

end;


end;
